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Abstract 

Much of the understanding of bulk liquids has progressed through study of the limiting case in 
which molecules interact via purely repulsive forces, such as a hard-core potential. In the same 
spirit, we report progress on the understanding of confined water by examining the behavior of 
water-like molecules interacting with planar walls via purely repulsive forces and compare our 
results with those obtained for Lennard-Jones (LJ) interactions between the molecules and the 
walls. Specifically, we perform molecular dynamics simulations of 512 water-like molecules which 
are confined between two smooth planar walls that are separated by 1.1 nm. At this separation, 
there are either two or three molecular layers of water, depending on density. We study two 
different forms of repulsive confinements, when the interaction potential between water-wall is 
(i) 1/r^ and (ii) WCA-like repulsive potential. We find that the thermodynamic, dynamic and 
structural properties of the liquid in purely repulsive confinements qualitatively match those for a 
system with a pure LJ attraction to the wall. In previous studies that include attractions, freezing 
into monolayer or trilayer ice was seen for this wall separation. Using the same separation as these 
previous studies, we find that the crystal state is not stable with 1/r^ repulsive walls but is stable 
with WCA-like repulsive confinement. However, by carefully adjusting the separation of the plates 
with 1 /r^ repulsive interactions so that the effective space available to the molecules is the same as 
that for LJ confinement, we find that the same crystal phases are stable. This result emphasizes 
the importance of comparing systems only using the same eff'ective confinement, which may differ 
from the geometric separation of the confining surfaces. 
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I. INTRODUCTION 



Confinement of water in nanopores affects many properties of water, such as freezing 



temperature, crystal structure , tlie glass transition temperature, and tlie 

position of tlie hypothesized liquid-liquid (LL) critical point 

mm 

■ Indeed, water 

confined in nanoscale geometries has received much recent attention, in part because of 
its importance in biology, engineering, geophysics and atmospheric sciences. The effects o_f 
different kinds of confinement have been studied, both using experiments and simulations 
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Bulk supercooled water - water cooled below the equilibrium freezing temperature 



shows many anomalous properties 



2G 
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221 . Experiments find that at low temperatures. 



various response functions, such as isothermal compressibility and specific heat, increase 
sharply. There has been comparatively less research on confined water. Using the ST2 



potential to model water confined between smooth plates 



231, a LL phase transition has been 



proposed. A liquid-to-amorphous transition is seen in simulations of water using the TIP4P 



potential 



24| confined in carbon nanotubes P]. Recent theoretical work 25| suggests that 



hydrophobic Lennard- Jones (LJ) confinement shifts the LL transition to lower temperature 
and lower pressure compared to bulk water, a feature also found in simulations of water 
confined between hydrophobic plates [iG^ . 

Confinement is known to enhance solidification of molecules that are more or less spher- 
ical H Q. Howeve. c„eM expen^en. o„ fl,.. of w.e. .Kow that wat. 
performs extremely well as a lubricant, suggesting that confined water may be more fluid 
than bulk water j2^. Recent experiments_show that water in hydrophilic confinement, when 

- a phenomenon also supported by simulation 
16l Il7l | show that hydrophobically con- 



cooled to v ery low T, does not freeze 
studies ISOlbll. In contrast, simulations 



fined water does freeze into different crystalline structures, which do not have counterparts 
in bulk water. Indeed monolayer, bilayer and trilayer ice have all been found in simula- 



tions 



17l |. Thus hydrophobic confinement seems to facilitate the freezing of 



water. However the reason for this facilitation is not yet fully understood. The hydrogen- 
bond interaction between water molecules is an order of magnitude stronger than Van der 
Waals attraction with the hydrophobic walls. Thus one may hypothesize that freezing in 
hydrophobic confinement depends critically on the separation between confining walls that 
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may distort or facilitate a particular crystalline structure, rather than on the weak details of 
the water-wall interaction potential. To test this hypothesis we perform molecular dynamics 
(MD) simulations of water in two different forms of repusive confinement. Specifically, we 
study: 

• The repulsive part of the LJ potential studied in Ref. 

• The same potential used in Ref. [3] but truncated and shifted such that there is no 
attractive part in the potential, analogous to the Weeks- Chandler- Anderson (WCA) 
potential [32]. This potential allows us to examine the role, if any, the attractive part 
of the water-wall LJ potential plays in determining the thermodynamics and structure 
of confined water. 

We compare the case when the water- wall interactions are purely repulsive ("repulsive 



confinement") with the studied case of pure LJ confinement jl6|. We also compare the 
freezing in repulsive confinements with the freezing found when the water-wall interactions 
are represented by an LJ interaction ("LJ hydrophobic confinement") 

This paper is organized as follows: In Sec. II, we provide details of our simulations and 
analysis methods. Simulation results for the liquid state are provided in Sec. III. In Sec. IV 
we discuss the freezing properties of our system. 

II. SIMULATION AND ANALYSIS METHODS 

We perform MD simulations of a system composed of water-like molecules confined be- 
tween two smooth walls. The molecules interact via the TIP5P pair potential 33i which, 
potential, treats each water molecule as a tetrahedral, rigid, and non- 
polarizable unit consisting of five point sites js^. The TIP5P potential predicts many of 
the anomalies of bulk water js^. For example, TIP5P reproduces the density anomaly at 
T = 277 K and P = 1 atm and its structural properties compare well with experiments 
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39|. TIP5P is known to crystallize at high pressures j36| within accessible 
computer simulation time scales, and shows a "nose-shaped" curve of temperature versus 
crystallization time 



a feature found in experimental data on water solutions . 
In our simulations, = 512 water molecules are confined between two smooth planar 
walls, as shown schematically in Fig. [T] The walls are located at 2;^ = ±0.55 nm (wall-wall 



separation of 1.1 nm), which results in ^ 2 — 3 layers of water molecules. Periodic boundary 
conditions are used in the x and y directions, parallel to the walls. 

We study two different forms of purely repulsive water-wall interaction. The first uses 
only the repulsive core, which we call the 1/r^ repulsive potential, 

" / \ 9' 

TT I \ A ' ^ 

U [z - zw) = 4eow 



\z — Zw\ 



(1) 



Here \z — zy^\ is the distance from the oxygen atom of a water molecule to the wall, while 
^ow = 0.25 kJ/mol and uow = 0.25 nm are potential parameters (Fig. |2I). Similar but 
different parameter values were used in previous confined water simulations using the TIP5P 
interaction potential Specifically, in Ref. the water- wall interaction was modelled 



using a 9-3 LJ potential with eow = 1-25 kJ/mol and gqw = 0.25 nm. We choose a different 
eow in the case of repulsive confinement so that the repulsion between the water and wall 
decays to almost zero where the 9-3 LJ -potential has a minimum. 

The second purely repulsive potential uses both attractive and repulsive terms of the 9-3 
L J potential, but truncates and shifts the potential at the position of the minimum to create 
a repulsive potential that exactly mimics the repulsion of ref. Q], in analogy to the WCA 
potential. 



U{z — Zw) = ' 



4eow[(^)^ - {^f] + if k - zw\ < 3^^ 



0"OW 



if |z — > 3^/^CroW5 



where eow = 1-25 kJ/mol and (Tqw = 0.25 nm. For the 1/r^ repulsive potential, we 
perform simulations for 56 state points, corresponding to seven temperatures T = 220 K, 
230 K, 240 K, 250 K, 260 K, 280 K, and 300 K, and eight "geometric densities" pg = 
0.60 g/cm^, 0.655 g/cm^, 0.709 g/cm^, 0.764 g/cm^, 0.818 g/cm^, 0.873 g/cm^, 0.927 g/cm^, 
and 0.981 g/cm^ - the same as studied in Ref. The geometric values of density do not 
take into account the fact that the repulsive interactions of molecules with the walls increases 
the overall amount of available space, since the eow parameter of the 1/r^ repulsive potential 
is smaller than the eow used for the LJ confined system |21 For systems confined by LJ 
interactions, there is a well-defined preferred distance from the wall, making it relatively 
straightforward to evaluate the "effective" density of molecules confined by the attractive 
wall. In our system with only repulsive interactions, there is no such preferred distance, as 
emphasized by Fig. |21 



We can approximate the effective density by examining the local density p{z) (Fig. 
We utilize the fact that p{z) has an inflection, and estimate the effective Lz by the loca- 
tion where the second derivative of p{z) = 0, or when first derivative of p{z) has a max- 
imum. We must also add to this value of the molecular diameter of water(0.278 nm) 
to calculate the real space available along 2;-direction. The resulting "effective densities" 
are p = 0.715 g/cm^, 0.777 g/cm^, 0.829 g/cm^, 0.890 g/cm^, 0.949 g/cm^, 1.000 g/cm^, 
1.060 g/cm^, and 1.115 g/cm^. We will use these effective densities throughout the paper, 
since they will be most comparable to the effective densities with LJ confinement. 

For the WCA potential, we perform simulations for 32 state points, correspond- 
ing to 8 different temperatures and 4 different "geometric densities", 0.60 g/cm^, 
0.655 g/cm^, 0.709 g/cm^, 0.764 g/cm^ respectively. These geometric densities correspond 
to "effective densities" 0.80 g/cm^, 0.88g/cm^, 0.95g/cm^, 1.02g/cm^ respectively. Note 
that these effective densities were calculated using the method described in 

We control the temperature using the Berendsen thermostat with a time constant of 5 
ps j4l| and use a simulation time step of 1 fs, just as in the bulk system Water- water 
interactions are truncated at a distance 0.9 nm as discussed in Ref. 0|. 



III. THERMODYNAMICS AND STRUCTURE 



One of the defining characteristics of water is the existence of a temperature of maximum 
density (TMD). Relative to bulk water - LJ confinement shifts the locus of the TMD to 
lower T by ~ 40 K Jjj. Additionally, the sharpness of the density maximum is markedly 
decreased in comparison to the bulk. Fig. |3] shows isochores of P for LJ confinement, 1/r^ 
repulsive confinement and WCA confinement, for similar densities. A TMD in this plot is 
coincident with the minimum in the isochore. For 1/r^ repulsive confinement, the minimum 
is very weak, but the location of the flatness in the isochore is near to that of the system with 
LJ confinement. This result suggests the the 1/r^ repulsive confinement further suppresses 
the structural ordering of the molecules that is known to be responsible to the presence of a 
density maximum. The TMD for the case of WCA confinement again appears at the same 
T as the 1/r^ repulsive confinement and LJ confinement cases but the isochore in the TMD 
region is flatter than for case the of LJ confinement confinemen. Hence both kinds of the 
repulsive confinement supprsess the structural ordering in lateral directions compared to the 
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case of bulk and LJ confinement. We further notice that the value of the lateral pressure 
P|l in the case of the L J confinement approaches the value of P\\ in case of WCA confinement 
at high temperatures. This behavior of P\\ for LJ confinement should be expected since at 
very high temperatures the molecules will not feel the potential minimum of the water-wall 
interaction. 

In order to compare the structural properties of repulsive confinement with those of 
LJ confinement, we calculate the lateral oxygen-oxygen radial distribution function (RDF) 
defined by 



(2) 



Here V is the volume, is the distance parallel to the walls between molecules i and j, Zi 
is the z-coordinate of the oxygen atom of molecule and 5 (a;) is the Dirac delta function. 
The Heaviside functions, Q[x\ restrict the sum to a pair of oxygen atoms of molecules 
located in the same slab of thickness bz = 0.1 nm. The physical interpretation of g\\{r) is 
that g\\{r)2iTrdr6z is proportional to the probability of finding an oxygen atom in a slab of 
thickness 6z at a distance r (parallel to the walls) from a randomly chosen oxygen atom. In 
a bulk liquid, this would be identical to g{r), the standard RDF. 

Figure IHl shows the temperature and density dependence of the lateral oxygen-oxygen pair 
correlation function for both repulsive (Fig. EK, Fig- Eb) and WCA (Fig. Efc, Fig. Ell) 
confinements. For both repulsive confinements, the qualitative behavior of the dependence 
of is the same. At low temperature and low density, the first two peaks in appear at 
r = 2.78 Aand r = 4.5 A, but at high densities the second peak moves to a larger distance. 
This behavior is nearly identical to that observed for water confined between LJ surfaces, 
and is discussed in detail in Ref. jlfi| . 

We also confirm the structural similarity with LJ confinement by calculating the lateral 
static structure factor S\\{q), defined as the Fourier transform of the lateral RDF g\\{r), 

^ll(^)-^E(e^'^-^^^~^^^^>- (3) 

Here the g-vector is the corresponding wave vector in the xy plane and r is the projection 
of the position vector on the xy plane. In Fig. El we show the temperature and pressure 
dependence of lateral structure factors for both repulsive confinements. For both forms of 
repulsive confinement, the temperature and density dependence of S\\ is similar. We find 



that confined water has a weaker pre-peak at ~ 18nm~^ compared to bulk water (Fig. [Tj), 
consistent with the possibility that the local tetrahedrality is weakened by repulsive confine- 
ments. Of the three forms of confinement, the for LJ confinement is most like bulk water 
(Fig. [7j). Local tetrahedrality becomes weaker in case of repulisive confinements compared 
to LJ confinement. Further, we see that the water in repulsive confinement is less 
structured in lateral directions compared to water in WCA confinement, indicated by a less 
sharp and broad pre-peak at ^ 18 nm^^ in S\i{q) (Fig. [7j). 



IV. FREEZING OF TIP5P WATER 



Bulk TIPS 
temperatures 



p 



water crystallizes within the simulation time for p ^ LIS g/cm^ at low 
Crystallization of confined water is seen in some simulations 0, Ig, Q|. 
A similar crystallization appears in simulations when an electric field is applied in lateral 
directions to a system of water confined between silica walls jisj]. 

At plate separation of 1.1 nm with hydrophobic LJ confinement, water crystallizes to 
trilayer ice ^|. From our simulations of TIPSP water in repulsive confinements with the 
same plate separation of 1.1 nm, we find that the system does not freeze within accessible 
simulation time scales for 1/r^ repulsive confinement; however the system freezes for WCA 
confinement. As a more stringent confirmation of this fact, we also use a starting ice config- 
uration obtained from simulations with LJ confinement for the same thickness, and confirm 
that the ice melts to a liquid with 1/r^ repulsive confinement. In Fig. |H1 we show the evo- 
lution of potential energy and lateral structure factor with time, when the crystal formed 
in LJ confinement Q| is kept between the 1/r^ repulsive walls. The potential energy first 
increases and then reaches its equilibrium value of the liquid potential energy accompanied 
by a structural change from a crystal (presence of sharp Bragg peaks) to a liquid (absence 
of Bragg peaks). 

Based on this observation, it is tempting to claim that repulsion inhibits crystallization, 
and that a preferable distance from the wall determined by the attractive portion of the LJ 
potential is necessary to induce crystallization. However, as discussed above for the same 
plate separation 1.1 nm, 1/r^ repulsive confinement with the chosen parameters increases 
the available space for molecules relative to LJ confinement. Hence to properly compare 
the crystallization behavior, we must adjust the separation of the wall so that the available 
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space for the water molecules is the same in both systems. We can make the available space 
the same by tuning the separation of the plates or by tuning the potential. By tuning the 
parameters (see Fig. |2l where the values of parameter eow and (Tqw for the modified 
repulsive potential are 1.25 kJ/mol and 0.23 nm respectively) of the 1 /r^repulsive potential 
in the repulsive system such that the density profile along the z-axis becomes similar, we 
have identical values of the available space between the plates (see Fig. El), and we find that 
an initial crystal configuration does not melt, emphasizing that the presence of the crystal 
is very sensitive to density and to plate separation - since the separation determines the 
accessible packing arrangements between the plates. Similar sensitivity to plate separation 



for monolayer ice was seen in ref. 15 1. 



In addition to examining the stability of initially crystalline structures, we also consider 
whether freezing from the liquid state occurs when we have the same effective plate separa- 
tion. We find that for the repulsively confined systems, the crystal will also spontaneously 
form if the available space between the plates is the same as that for which initially crys- 
talline configurations are stable. Hence plate separation appears to be the dominant cause 
in determining whether or not a crystal will form. 
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FIG. 1: Perspective view of the system, showing the 512 water molecules confined between two walls 
perpendicular to the z-direction. Note that the confining plates are located along the z-direction 
and are separated by 2-3 molecular layers of water. 
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FIG. 2: Water-wall interaction potentials namely the 9-3 LJ, repulsive, and WCA potential 
as a function of distance of water molecules \z — Zu]\ from the center of one of the walls (shaded 
rectangle) . 
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FIG. 3: Density profile p{z) along z-direction for four different bulk densities at T=250 K for 
repulsive confinement, (b) Density profile p{z) along z-direction for four different bulk densities at 



T=220 K for the case of WCA confinement. Both the lepu. 



sive confinements show similar layering 



of water molecules as seen in the case of L J confinement l£ , ll 
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FIG. 4: Lateral pressure Py for one isochore for the purely repulsive, LJ confinement and WCA 
cases. Here the effective density is p = 0.829 g/cm^ for the repulsive potential and p = 
0.950 g/cm'^ for LJ confinement. These effective densities for both systems correspond to the same 
geometric density of pg = 0.709 g/cm'^. All forms of confinement show a TMD, indicated by the 
minimum of the pressure; however the TMD is very "flat" for repulsive and WCA confinement. 
As expected the value of Py approaches the value of Py for the case of LJ confinement at high 
temperatures. 
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FIG. 5: Lateral oxygen-oxygen pair correlation function g\\{r) for the case of repusive and 
WCA confinements. Shown are four different densities at two fixed temperatures (a) T = 220 K 
and (b) T = 300 K for repulsive confinement and (c) and (d) for WCA confinement. Note 
that with increasing density, the second neighbor peak at ~ 0.45 nm becomes less pronounced and 
at high temperatures moves to a larger distance seen in LJ confinement [l^ . 
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FIG. 6: Lateral structure factor for the case of repulsive and WCA confinements. Shown 
are four different densities at two fixed temperatures (a) T = 220 K and (b) T = 300 K for 
repulsive confinement and (c) and (d) for WCA confinement for the same temperatures. The first 
peak of 5|| {q) corresponding to the hydrogen-bonds weakens as density is increased and is absent 
at high densities and high temperatures. The first peak of S^{q) in the case of repulsive 
confinement weaker than the LJ confinement. 
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FIG. 7: A comparison of the structure factors for different confinements with bulk water structure 
factor at T=250K. For the comparion we choose the effective densities in confinements close to each 
other. Densities for repulsive confinement, LJ confinement, WCA-confinement are chosen 
to be 0.950 g/cm^. We choose the structure factor for bulk water at density 1.00 g/cm'^ . A 
diminished peak at ~ 18 nm^^ shows that the local tetrahedral structure is weakened in case of all 
forms of confinement. A further comparison of LJ confinement at 0.950 g/cm^ with 1/r^ repulsive 
confinement shows that water in 1/r^ repulsive confinement is less tetrahedral, as the first peak of 
S^^{q) is much weaker than the first peak of S\^{q) for LJ confinement. Water in WCA confinement 
is more structured than 1/r^ repulsive confinement, but is less structured than LJ confinement. 
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FIG. 8: (a) A plot of potential energy as a function of t, when the crystal formed in LJ confine- 
ment 16] is kept between the repulsive wall, (b) The crystal structure indicated by the sharp 
Bragg peaks melts, and (c) turns into a liquid indicated by the absence of Bragg peaks. 
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FIG. 9: Density profile p{z) of water along z-direction for different potentials at the same geometric 
density. The repulsive confinement system freezes spontaneously when the parameters of the 
potential are modified such that the effective Lz calculated from the p{z) (red dotted line) for the 
repulsive system is same as that for the LJ system (blue-dashed line) (see Fig. [21). 
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FIG. 10: (a) Potential energy as a function of time t, for 1/r^ repulsive confinement when the 
effective Lz is same as the effective Lz of LJ confinement at T = 260 K and geometric density 
Pg = 0.981 g/cm^. The confined water spontaneously freezes, indicated by the drop in potential 
energy, (b) The structure factor of the ice such formed resembles the trilayer ice seen in case of LJ 
confinement ^, 1^ . 
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